
function [I,P1] = ewrlnc_2layers(nb_sims,I_block,P1_block,m,nb_transmissions,R1,R2)

% Description:
% This function simulates Expanding Windows Random Linear Network Coding
%
% Input Parameters:
% m                     is the field parameter so thath 2^m = elements in field
% nb_transmissions      is the number of transmitted packets (20)
% nb_sims               is the number of sims to run / received pkt's
% I_block               is the number of symbols in layer 1
% P1_block              is the number of symbols in layer 2
% P2_block              is the number of symbols in layer 3
% R1                    is the probability of sending layer 1
% R2                    is the probability of sending layer 2
% R3                    is the probability of sending layer 3
%
% Output Parameters:
% I                     is the number of times I block was decoded
% P1                    is the number of times P1 block was decoded
% P2                    is the number of times P1 block was decoded
%
% To do:
%       - Nothing?
%
% Profiler:
%       - rank(gf(... )) takes all the time.

%% Setup

R=[R1 R2]; % The discrete layer probability distribution is given as
I=zeros(nb_sims,1); % nb of I recv, malloc
P1=zeros(nb_sims,1); % nb of P1 recv, malloc
% P2=zeros(nb_sims,1); % nb of P2 recv, malloc

% Dynamic setup of the three layers based on input parameters.
I_pkts=I_block; % nb of packets in I block
P1_pkts=P1_block; % nb of packets in P1 block
% P2_pkts=P2_block; % nb of packets in P2 block
I_P1_pkts=I_pkts+P1_pkts; % nb of packets in I + P1 block
% I_P1_P2_pkts=I_pkts+P1_pkts+P2_pkts; % nb of packets in I + P1 + P2 block
gen_size=I_pkts+P1_pkts; % nb of packets in generation = generation size.

%% Simulation

for sim=1:nb_sims % Repeat for nb_sims
    
    % Allocate M_recv matrix
    M_recv_1=zeros(nb_transmissions,gen_size); % Initial M_recv with layer 1 entries
    M_recv_2=zeros(nb_transmissions,gen_size); % Initial M_recv with layer 2 entries
%     M_recv_3=zeros(nb_transmissions,gen_size); % Initial M_recv with layer 3 entries
    
    % layer counters, index into M_recv matrices
    layer_1_index=0;
    layer_2_index=0;
%     layer_3_index=0;
    
    for j=1:nb_transmissions % Repeat for number of transmissions/gen
        
        r=rand(1); % Decide on this random draw what layer to send
        
        if r<R(1)
            % Layer 1 should be sent
            layer_1_index=layer_1_index+1;
            new_row=[randi([0 2^m-1],1,I_pkts) zeros(1,P1_pkts)];
            M_recv_1(layer_1_index,:)=new_row;
        else
            % Layer 2 should be sent
            layer_2_index=layer_2_index+1;
            new_row=[randi([0 2^m-1],1,I_P1_pkts)];
            M_recv_2(layer_2_index,:)=new_row;
        end
        
    end
    
    % Collect stats on decoded data packets
    % Assumption: rank(gf(Matrix)) can be trusted!!!!!
    
    % Speed: rank(gf(... )) takes 3/4 of the time spent in this script
    % Speed: Short circuit the rank(gf(...)) depending on nb_transmissions
    
    % NOTE: We should count layers, not block of data?
    %       commenting out in the following!....
    
    %     if nb_transmissions>=gen_size && rank(gf([M_recv_1;M_recv_2;M_recv_3],m))>=gen_size
    %         % I,P1 and P2 can be decoded
    % %         I(sim)=I(sim)+1;
    % %         P1(sim)=P1(sim)+1;
    %         P2(sim)=P2(sim)+1;
    %     elseif nb_transmissions>=I_P1_pkts && rank(gf([M_recv_1;M_recv_2],m))>=I_P1_pkts
    %         % I and P1 can be decoded
    % %         I(sim)=I(sim)+1;
    %         P1(sim)=P1(sim)+1;
    %     elseif nb_transmissions>=I_pkts && rank(gf(M_recv_1,m))>=I_pkts
    %         % I decoded can be decoded
    %         I(sim)=I(sim)+1;
    %     end
    %
    
    % The >new< way of checking rank.... (BK,EV)
    
    if nb_transmissions>=I_P1_pkts && rank(gf([M_recv_1;M_recv_2],m))>=I_P1_pkts
        % I and P1 can be decoded
        %         I(sim)=I(sim)+1;
        P1(sim)=P1(sim)+1;
    end
    
    if nb_transmissions>=I_pkts && rank(gf(M_recv_1,m))>=I_pkts
        % I decoded can be decoded
        I(sim)=I(sim)+1;
    end
    
    
end

I=sum(I); % How many times was I block decoded
P1=sum(P1); % How many times was P1 block decoded
% P2=sum(P2); % How many times was P2 block decoded

% Returning ...

end





